/* ThermalRegulation�Ĳ���time�ɴ����ø�ΪΪ��ֵ
 * 
 */
package com.simulation;

class T25NodeBody
{
    // the controlled system
    double m_heatCapacity[] = new double[25];  //C the heat capacitance of each part of body
    double m_tBody[] = new double[25]; // T the temperature of body
    double m_tSkin[] = new double[6]; // the temperature of skin
    double m_tCore[] = new double[7];
    double m_tChangeRate[] = new double[25]; //F rate of tempereture change in compart n
    double m_heatFlowRate[] = new double[25]; // HF rate of heat flow into or from n
    double m_thermoConductivity[] = new double[24]; //TC thermal conductivity between part n and n+1
    double m_heatTransferRate[] = new double[24]; //TD convective heat transfer rate n to n+1
    double m_pSkin[] = new double[6]; //the pressure of the skin
    double m_cWaterVapor[] = new double[6];

    double m_bHeatMetabolic[] = new double[24];  //QB the basic heat produced by metabolic
    double m_bHeatEvaporation[] = new double[24]; //EB the basic heat loss from evaporation
    double m_heatConvectionCenter[] = new double[24]; //BC convective heat transfer between central blood and n
    double m_bBloodFlow[] = new double[24];    //BFB basic effective blood flow to
    double m_totalHeatMetabolic[] = new double[24]; // Q total metabolic heat production in n
    double m_totalHeatEvaporation[] = new double[24]; // E total evaporative heat loss from n
    double m_totalBloodFlow[] = new double[24];       // BF total effective blood flow to n
    double m_heatConvectionCoe[] = new double[6]; //HC the environmental convection heat transfer coefficient for segement i
    double m_surface[] = new double[6]; //S the surface of segement i
    double m_heatTransferCoe[] = new double[6]; //H total environmental heat transfer coefficient for segment i
    double m_vaporPressureEnvironment[] = new double[6];// PCL
    double m_vaporPressureSkin[] = new double[6];       // PSKIN
    double m_vaporPressureSatuated[] = new double[6];   // PSAT
    double m_heatResSkin; //RES the heat resistance of evaporation of skin
    double m_heatResAir[] = new double[6];//REA the heat resistance of evaporation in the air layer

    double m_sweatRate[] = new double[6]; //CMSKIN sweating rate
    double m_sweatingRateSkin[] = new double[6]; // MRSW
    double m_mAccumulatedSkin[] = new double[6]; //      the moisture volume accumulated in the skin surface
    double m_sAccumulatedSkin[] = new double[6]; //HH    the sweat volume accumulated in the skin surface
    double m_eSkin[] = new double[6]; // Esk, the heat loss of evaporation of skin

    // control coefficients of contrlling system
    double m_sweatHeadCore;//CSW sweating from head core
    double m_sweatSkin; //SSW sweating from skin
    double m_sweatSkinAndHead; // PSW
    double m_vasodHeadCore; //CDIL vasodilatation from skin
    double m_vasodSkin; // SDIL vasodilatation from skin;
    double m_vasodSkinAndHead; //PDIL
    double m_vasocHeadCore;// CCON vasoconstriction from head core;
    double m_vasocSkin; // SCON vasoconstriction from skin
    double m_vasocSkinAndhead;// PCON
    double m_shiveringHead; //CCHIL shivering from head core
    double m_shiveringSkin; // SCHIL shivering from skin
    double m_shiveringSkinAndHead;// PCHIL

    // used in the description of the controlling system
    double m_temperatureReference[] = new double[25];//TSET set point or reference point for receptors in compart i
    double m_outputError[] = new double[25]; // ERROR output from thermoreceptors in compartment i
    double m_dynSensitivity[] = new double[25]; //RATE  dynamic sensitivity of thermareceptors in compart n
    double m_outputWarm[] = new double[25]; //WORM output form warm receptors in n
    double m_outputCold[] = new double[25]; //COLD output from cold receptors in n
    double m_integratedOutputWarm;  //WORMS  integrated output from warm receptors
    double m_integratedOutputCold;  //COLDS  integrated output from cold receptors
    double m_totalSweat; //SWEAT  total efferent sweat command
    double m_totalVascodSkin; // DILAT total efferent skin vasodilatation command
    double m_totalVasocSkin;  // STRIC  total efferent skin vasoconstriction command
    double m_totalShivering;  // CHILL  total efferent shivering command

    double m_receptorSkin[] = new double[6]; // SKINR fraction of skin receptors in segment i
    double m_sweatCommandSkin[] = new double[6]; //SKINS fraction of sweating command applicable to skin of segment i
    double m_vasodCommandSkin[] = new double[6]; //SKINV fraction of vasodilatation commmand applicable to skin of segment i
    double m_vasocCommandSkin[] = new double[6]; //SKINC fraction of vasoconstriction command applicable to skin of segment i
    double m_muscleWork[] = new double[6]; //MWORK fraction of total work done by muscles in segnment i
    double m_muscleShivering[] = new double[6]; //MCHIL fraction of total shivering occurring in muscles in segnment i
    double m_work;
    double m_rh[] = new double[6];
	
	T25NodeBody() {
		for(int i=0;i<24;i++)
		{
			m_totalHeatEvaporation[i] = 0;
		}
	};
	
    void InitialValue(BodyType body) {
    	for(int i=0;i<25;i++)
        {
            m_heatCapacity[i]=body.C[i];// C[i];
            m_temperatureReference[i]=body.TSET[i];// TSET[i];
            m_tBody[i]=body.T[i];// T[i];
        }
        for(int i=0;i<24;i++)
        {
            m_bHeatMetabolic[i]=body.QB[i];
            m_bHeatEvaporation[i]=body.EB[i];
            m_bBloodFlow[i]=body.BFB[i];
            m_thermoConductivity[i]=body.TC[i];
        }
        for(int i=0;i<6;i++)
        {
            m_surface[i]=body.S[i];
            m_heatConvectionCoe[i]=body.HC[i];

            m_receptorSkin[i]=body.SKINR[i] ;
            m_sweatCommandSkin[i]=body.SKINS[i] ;
            m_vasodCommandSkin[i]=body.SKINV[i];
            m_vasocCommandSkin[i]=body.SKINC[i] ;
            m_muscleWork[i]=body.MWORK[i] ;
            m_muscleShivering[i]=body.MCHILL[i];
        }

        m_sweatHeadCore=body.CSW;// 372.0; //CSW
        m_sweatSkin=body.SSW;// 33.7; //SSW
        m_sweatSkinAndHead=body.PSW;// 0.0; //PSW
        m_vasodHeadCore=body.CDIL;// 136.0; //CDIL
        m_vasodSkin=body.SDIL;// 17.0;  //SDIL
        m_vasodSkinAndHead=body.PDIL;// 0.0; //PDIL
        m_vasocHeadCore=body.CCON;// 10.8;  //CCON
        m_vasocSkin=body.SCON;// 10.8;  //SCON
        m_vasocSkinAndhead=body.PCON;// 0.0; //PCON
        m_shiveringHead=body.PCON;// 13.0; //CCHIL
        m_shiveringSkin=body.SCHIL;// 0.40; //SCHIL
        m_shiveringSkinAndHead=body.PCHIL;// 1.0;  //PCHIL

        m_heatResSkin=body.Res;// 333.3;  //RES

        for(int i=0;i<25;i++)
        {
            m_dynSensitivity[i]=0.0;
            m_tChangeRate[i]=0.0;
        }

        m_work=0.0;

        //DELT=1;

        if(m_work-86.5<=0)  // the value of work must be more than 86.5
        {
            m_work=0.0;
        }
        else
        {
            m_work=m_work-86.5;
        }

        for(int i=0;i<6;i++)
        {
            m_sweatRate[i]=0.0;
            //SKIN_T[i] = 33.0+273.16;
            m_tSkin[i]=m_tBody[4*i+3]+273.16; //33.0+273.16;
            m_cWaterVapor[i] =TPhysicalPro.WaterVaporConcentration(m_tSkin[i])*0.6;
            m_eSkin[i] = 0.0;
            m_mAccumulatedSkin[i]=0;
            m_rh[i]=0.6;
        }
    };
    
    void ThermalRegulation(double tClothing[], double rhClothing[],double tEnvironment,double rhEnvironment,double time, int itime,double iRatio[]) {
    	double pClothing[] = new double[6],pEnvironment;
    	double sweat,dilat,strict,chill;
    	double pSatuated[] = new double[6];
    	double totalHeatMetabolic;
    	double dt;
    	int LTIME;

    	do
    	{
    		pClothing[0] = rhClothing[0] * Math.exp(23.48 - 3996.0 / (tClothing[0] - 39.1));
    		pClothing[1] = rhClothing[1] * Math.exp(23.48-3996.0/(tClothing[1]-39.1));
    		pClothing[2] = rhClothing[2] * Math.exp(23.48 - 3996.0 / (tClothing[2] - 39.1));
    		pClothing[3] = rhClothing[3] * Math.exp(23.48 - 3996.0 / (tClothing[3] - 39.1));
    		pClothing[4] = rhClothing[4] * Math.exp(23.48 - 3996.0 / (tClothing[4] - 39.1));
    		pClothing[5] = rhClothing[5] * Math.exp(23.48 - 3996.0 / (tClothing[5] - 39.1));

    		pEnvironment=rhEnvironment*Math.exp(23.48-3996.0/(tEnvironment-39.1));

    		for(int i=0;i<25;i++)//使用25节点模型
    		{
    			m_outputWarm[i]=0.0;
    			m_outputCold[i]=0.0;
    			m_outputError[i]=m_tBody[i]-m_temperatureReference[i]+m_dynSensitivity[i]*m_tChangeRate[i];

    			if(m_outputError[i]<0.0)
    			{
    				m_outputCold[i]=-m_outputError[i];
    			}
    			else if(m_outputError[i]>0.0)
    			{
    				m_outputWarm[i]=m_outputError[i];
    			}
    		}

    		m_integratedOutputWarm=0.0;
    		m_integratedOutputCold=0.0;
    		for(int i=0;i<6;i++)
    		{
    			m_integratedOutputWarm=m_integratedOutputWarm + m_outputWarm[4*i+3]*m_receptorSkin[i];
    			m_integratedOutputCold=m_integratedOutputCold + m_outputCold[4*i+3]*m_receptorSkin[i];
    		}

    		m_totalSweat=m_sweatHeadCore*m_outputError[0]+2.0*m_sweatSkin*(m_integratedOutputWarm-m_integratedOutputCold)
    			+m_sweatSkinAndHead*m_outputError[0]*(m_integratedOutputWarm-m_integratedOutputCold);
    		m_totalVascodSkin=m_vasodHeadCore*m_outputError[0]+m_vasodSkin*(m_integratedOutputWarm-m_integratedOutputCold)
    			+m_vasodSkinAndHead*m_outputWarm[0]*m_integratedOutputWarm;
    		m_totalVasocSkin=-m_vasocHeadCore*m_outputError[0]-m_vasocSkin*(m_integratedOutputWarm-m_integratedOutputCold)
    			+ m_vasocSkinAndhead*m_outputCold[0]*m_integratedOutputCold;
    		m_totalShivering=(m_shiveringHead*m_outputError[0]+m_shiveringSkin*(m_integratedOutputWarm-m_integratedOutputCold))
    			*m_shiveringSkinAndHead*(m_integratedOutputWarm-m_integratedOutputCold);

    		if(m_totalSweat<=0.0)
    		{
    			m_totalSweat=0.0;
    		}
    		if(m_totalVascodSkin<=0.0)
    		{
    			m_totalVascodSkin=0.0;
    		}
    		if(m_totalVasocSkin<=0.0)
    		{
    			m_totalVasocSkin=0.0;
    		}
    		if(m_totalShivering<=0.0)
    		{
    			m_totalShivering=0.0;
    		}

    		for(int i=0;i<6;i++)
    		{
    			m_sweatingRateSkin[i]=m_sweatCommandSkin[i]*m_totalSweat/0.7;

    			pSatuated[i]=Math.exp(23.48-3996.0/(m_tBody[4*i+3]+273.16-39.1));

    			//重要，更新m_pSkin，这个参数在计算最终湿度时会用到
    			if(Math.abs(m_mAccumulatedSkin[i])<1.0e-6)
    			{
    				m_pSkin[i]=(pSatuated[i]*m_heatResAir[i]+pClothing[i]*m_heatResSkin
    					+m_sweatingRateSkin[i]*0.7*m_heatResAir[i]
    				*m_heatResSkin/m_surface[i])/(m_heatResAir[i]+m_heatResSkin);

    				if((m_pSkin[i]>pSatuated[i])||(m_pSkin[i]<0))
    				{
    					m_pSkin[i]=pSatuated[i];
    				}
    			}
    			else
    			{
    				m_pSkin[i]=pSatuated[i];

    				if(m_mAccumulatedSkin[i]>35.0*m_surface[i])
    				{
    					m_mAccumulatedSkin[i]=35.0*m_surface[i];
    				}
    			}

    			m_sAccumulatedSkin[i]=m_sweatingRateSkin[i]+(pSatuated[i]-m_pSkin[i])*m_surface[i]/(m_heatResSkin*0.7)
    				-(m_pSkin[i]-pClothing[i])*m_surface[i]/(m_heatResAir[i]*0.7);

    		}
    		
    		for(int i=0;i<6;i++)//计算热代谢、血流量、热蒸发
    		{
    			//内核
    			m_totalHeatMetabolic[4*i]=m_bHeatMetabolic[4*i];
    			m_totalBloodFlow[4*i]=m_bBloodFlow[4*i];
    			m_totalHeatEvaporation[4*i]=0.0;

    			//肌肉
    			m_totalHeatMetabolic[4*i+1]=m_bHeatMetabolic[4*i+1]+m_muscleWork[i]*m_work+m_muscleShivering[i]*m_totalShivering;
    			m_totalBloodFlow[4*i+1]=m_bBloodFlow[4*i+1]+m_totalHeatMetabolic[4*i+1]-m_bHeatMetabolic[4*i+1];
    			m_totalHeatEvaporation[4*i+1]=0.0;

    			//脂肪
    			m_totalHeatMetabolic[4*i+2]=m_bHeatMetabolic[4*i+2];
    			m_totalBloodFlow[4*i+2]=m_bBloodFlow[4*i+2];
    			m_totalHeatEvaporation[4*i+2]=0.0;

    			//皮肤
    			double ratio = iRatio[i]/100;//这里用到了衣物覆盖率
    			m_totalHeatMetabolic[4*i+3]=m_bHeatMetabolic[4*i+3];
    			m_totalBloodFlow[4*i+3]=(m_bBloodFlow[4*i+3]+m_vasodCommandSkin[i]*m_totalVascodSkin)/(1.0+m_vasocCommandSkin[i]*m_totalVasocSkin);
    			m_totalHeatEvaporation[4*i+3]=(m_pSkin[i]-pClothing[i])*m_surface[i]*ratio/m_heatResAir[i] + (1-ratio)*(m_pSkin[i]-pSatuated[i]*rhClothing[i])*m_surface[i]/m_heatResAir[i];
    	
    			m_eSkin[i] = m_totalHeatEvaporation[4*i+3]/m_surface[i];

    		}

    		totalHeatMetabolic=0.0;

    		for(int i=0;i<24;i++)
    		{
    			totalHeatMetabolic=totalHeatMetabolic+m_totalHeatMetabolic[i];
    		}
    		//这里用到外界温度
    		m_totalHeatEvaporation[4]=(0.0014*(34.0-(tEnvironment-273.16))+0.017*(5.867-pEnvironment/1000.0))*totalHeatMetabolic;

    		for(int i=0;i<24;i++)
    		{
    			m_heatConvectionCenter[i]=m_totalBloodFlow[i]*(m_tBody[i]-m_tBody[24]);
    			m_heatTransferRate[i]=m_thermoConductivity[i]*(m_tBody[i]-m_tBody[i+1]);
    		}
    		
    		for(int i=0;i<6;i++)
    		{
    			m_heatFlowRate[4*i]=m_totalHeatMetabolic[4*i]-m_totalHeatEvaporation[4*i]-m_heatConvectionCenter[4*i]-m_heatTransferRate[4*i];
    			m_heatFlowRate[4*i+1]=m_totalHeatMetabolic[4*i+1]-m_heatConvectionCenter[4*i+1]+m_heatTransferRate[4*i]-m_heatTransferRate[4*i+1];
    			m_heatFlowRate[4*i+2]=m_totalHeatMetabolic[4*i+2]-m_heatConvectionCenter[4*i+2]+m_heatTransferRate[4*i+1]-m_heatTransferRate[4*i+2];
    			m_heatFlowRate[4*i+3]=m_totalHeatMetabolic[4*i+3]-m_heatConvectionCenter[4*i+3]-m_totalHeatEvaporation[4*i+3]+m_heatTransferRate[4*i+3]-m_heatTransferCoe[i]*(m_tBody[4*i+3]-(tClothing[i]-273.16));
    		}

    		m_heatFlowRate[24]=0.0;
    		for(int i=0;i<24;i++)
    		{
    			m_heatFlowRate[24]=m_heatFlowRate[24]+m_heatConvectionCenter[i];
    		}

    		m_heatFlowRate[24]=m_heatFlowRate[24]-0.08*m_work;

    		//设置仿真步长
    		dt = 1/360.0;
    		for(int i=0;i<24;i++)
    		{
    			m_tChangeRate[i]=m_heatFlowRate[i]/m_heatCapacity[i];

    			if((Math.abs(m_tChangeRate[i])*dt-0.01)>0.0)
    			{
    				dt=0.01/Math.abs(m_tChangeRate[i]);
    			}
    		}
    		m_tChangeRate[24] = m_heatFlowRate[24] / m_heatCapacity[24];
    		
    		//重要，计算m_tBody，这个参数在计算最终温度时会用到。前面计算m_tChangeRate、m_heatFlowRate等，其实都是为了求出m_tBody.
    		for(int i=0;i<24;i++)
    		{
    			m_tBody[i]=m_tBody[i]+m_tChangeRate[i]*dt;
    			if (i%4==3)
    			{
    				if (m_tBody[i]<30)
    				{
    					m_tBody[i]=30;
    				}
    			}
    		}
    		m_tBody[24] = m_tBody[24] + m_tChangeRate[24] * dt;
    		
    		for(int i=0;i<6;i++)
    		{
    			m_mAccumulatedSkin[i]=m_mAccumulatedSkin[i]+m_sAccumulatedSkin[i]*dt;
    			if(m_mAccumulatedSkin[i]<0.0)
    			{
    				m_mAccumulatedSkin[i]=0.0;
    			}
    		}
    		time = time+dt;   
    		LTIME = (int)(Math.floor(360.0*time));
    	}
    	while((LTIME-1-itime)<0);

    	for(int i=0;i<6;i++)
    	{
    		m_tSkin[i]=m_tBody[4*i+3]+273.16;//皮肤温度
    		m_tCore[i]=m_tBody[4*i];//内核温度
    		m_cWaterVapor[i] =m_pSkin[i]*0.018/(8.31*m_tSkin[i])*1.0e-3;//水蒸汽
    		m_rh[i]=m_pSkin[i]/TPhysicalPro.SaturatedWaterVaporP(m_tSkin[i]);//湿度
    	}
    };

};
